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METHOD FOR DETERMINING PARAMETERS OF EARTH 
FORMATIONS SURROUNDING A WELL BORE 

BACKGROUND OF THE INVENTION 

5 Technical Field of the Invention 

The present invention relates to well logging, and more particularly, to 
a method for determining formation parameters around a well bore. 

Description of Related Art 

Modern petroleum drilling and production operations demand a great 

10 quantity of information relating to parameters and conditions downhole. Such 
information typically includes characteristics of the earth formations traversed 
by the wellbore, in addition to data relating to the size and configuration of the 
borehole itself Oil well logging has been known in the industry for many 
years as a technique for providing information to a formation evaluation 

15 professional or driller regarding the particular earth formation being drilled. 
The collection of information relating to conditions downhole, which commonly 
is referred to as "logging," can be performed by several methods. These 
methods include measurement while drilling, MWD, and logging while drilling, 
LWD, in which a logging tool is carried on a drill string during the drilling 

20 process. The methods also include wireline logging. 

in conventional oil well wireline logging, a probe or "sonde" is lowered 
into the borehole after some or all of the well has been drilled, and is used to 
determine certain characteristics of the formations traversed by the borehole. 
The sonde may include one or more sensors to measure parameters 

25 downhole and typically is constructed as a hermetically sealed cylinder for 
housing the sensors, which hangs at the end of a long cable or "wireline." 
The cable or wireline provides mechanical support to the sonde and also 
provides electrical connections between the sensors and associated 
instrumentation within the sonde, and electrical equipment located at the 

30 surface of the well. Normally, the cable supplies operating power to the 
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sonde and is used as an electrical conductor to transmit information signals 
from the sonde to the surface. In accordance with conventional techniques, 
various parameters of the earth's formations are measured and correlated 
with the position of the sonde in the borehole as the sonde is pulled uphole. 

5 A chart or plot of an earth parameter or of a logging tool signal versus 

the position or depth in the borehole is called a "log." The depth may be the 
distance from the surface of the earth to the location of the tool in the 
borehole or may be true depth, which is the same only for a perfectly vertical 
straight borehole. The log of the tool signal or raw data often does not 

10 provide a clear representation of the earth parameter which the formation 
evaluation professional or driller needs to know. The tool signal must usually 
be processed to produce a log which more clearly represents a desired 
parameter. The log is normally first created in digital form by a computer and 
stored in computer memory, on tape, disk, etc. and may be displayed on a 

1 5 computer screen or printed in hard copy form. 

The sensors used in a wireline sonde usually include a source device 
for transmitting energy into the formation, and one or more receivers for 
detecting the energy reflected from the formation. Various sensors have been 
used to determine particular characteristics of the formation, including nuclear 

20 sensors, acoustic sensors, and electrical sensors. See generally J. Lab, A 
Practical Introduction to Borehole Geophysics (Society of Exploration 
Geophysicists 1986); D.R. Skinner, Introduction to Petroleum Production, 
Volume 1, at 54-63 (Gulf Publishing Co. 1981). 

For a formation to contain petroleum, and for the formation to permit the 

25 petroleum to flow through it, the rock comprising the formation must have 
certain well-known physical characteristics. One characteristic is that the 
formation has a certain range of measurable resistivity (or conductivity), which 
in many cases can be determined by inducing an alternating electromagnetic 
field into the formation by a transmitter coil arrangement. The 

30 electromagnetic field induces alternating electric (or eddy) currents in the 
formation in paths that are substantially coaxial with the transmitter. These 
currents in turn create a secondary electromagnetic field in the medium, 
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inducing an alternating voltage at the receiver coil. If the current in the 
transmitter coil is kept constant, the eddy current intensity is generally 
proportional to the conductivity of the formation. Consequently, the 
conductivity of the formation determines the intensity of the secondary 
5 electromagnetic field, and thus, the amplitude of the voltage at the receiver 
coil. See generally, James R. Jordan, et al., Well Logging II -- Electric And 
Acoustic Logging, SPE Monograph Series, Volume 10, at 71-87 (1986). 

An exemplary induction tool is shown in the prior art drawing of Figure 
1 , in which one or more transmitters (T) and a plurality of receivers (Ri) are 

10 shown in a logging sonde. Each transmitter or receiver may be a set of coils, 
with modern array induction tools having several receivers, e.g. R1, R2, R3, 
and R4, of increasing transmitter-to-receiver spacing to measure 
progressively deeper into the formation. 

In a conventional induction tool such as that shown in Figure 1, the coils 

15 are wound coaxially around a cylindrical mandrel. Both transmitter coils and 
receiver coils are solenoidal, and are wound coaxial with the mandrel. Such 
coils would therefore be aligned with the principal axis of the logging tool, 
which is normally also the central axis of the borehole and is usually referred 
to as the z-axis. That is, the magnetic moments of the coils are aligned with 

20 the axis of the mandrel on which they are wound. The number, position, and 
numbers of turns of the coils are arranged to null the signal in a vacuum due 
to the mutual inductance of transmitters and receivers. 

During operation, an oscillator supplies alternating current to the 
transmitter coil or coils, thereby inducing current in the receiver coil or coils. 

25 The voltage of the current induced in the receiver coils results from the sum of 
all eddy currents induced in the surrounding formations by the transmitter 
coils. Phase sensitive electronics measure the receiver voltage that is in- 
phase with the transmitter current divided by magnitude of the transmitter 
current When normalized with the proper scale factor, this provides signals 

30 representing the apparent conductivity of that part of the formation through 
which the transmitted signal passed. The out-of-phase, or quadrature, 
component can also be useful because of its sensitivity to skin effect although 
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it is less stable and is adversely affected by contrasts in the magnetic 
permeability. 

As noted, the induced eddy currents tend to flow in circular paths that 
are coaxial with the transmitter coil. As shown in Figure 1, for a vertical 
5 borehole traversing horizontal formations, there is a general symmetry for the 
induced current around the logging tool. In this ideal situation, each line of 
current flow remains in the same formation along its entire flow path, and 
never crosses a bed boundary. 

In many situations, as shown for example in Figure 2, the wellbore is 

10 not vertical and the bed boundaries are not horizontal. The well bore in 
Figure 2 is shown with an inclination angle 0 measured relative to true 
vertical. A bed boundary between formations is shown with a dip angle a. 
The inclined wellbore strikes the dipping bed at an angle (J. As a result, the 
induced eddy currents flow through more than one media, encountering 

15 formations with different resistive properties. The resulting logs are distorted, 
especially as the dip angle a of the bed boundaries increases. If the logging 
tool traverses a thin bed, the problem becomes even more exaggerated. 

As shown in the graph of Figure 3A, an induction sonde traversing a 
dipping bed produces a log with distortions normally referred to as "horns". 

20 The more severe the dip angle, the less accurate is the measurement with 
depth. Figure 3A represents a computer simulation of a log that would be 
generated during logging of a ten-foot thick bed (in actual depth), with 
different plots for different dip angles. Figure 3B shows a computer simulation 
of a log which would be generated if the thickness of the bed were true 

25 vertical depth, with different plots for different dip angles. As is evident from 
these simulated logs, as the dip angle increases, the accuracy and 
meaningfulness of the log decreases. In instances of high dip angles, the 
plots become virtually meaningless in the vicinity of the bed boundaries. 

Figures 3A and 3B also illustrate that even for a vertical well traversing 

30 horizontal formations, the actual electrical signal or data produced by an 
induction logging tool is quite different from an exact plot of formation 
resistivities. In these figures the desired representations of formation 
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resistivity are the dashed line square wave shapes 10 and 20. The actual 
resistivity within a layer is generally uniform so that there are abrupt changes 
in resistivity at the interfaces between layers. However, logging toots have 
limited resolution and do not directly measure these abrupt changes. When 
5 the transmitter coil T in Figure 1 is near an interface, as illustrated, its 
transmitted signal is split between layers of differing resistivity. As a result, 
the raw data or signal from the logging tool is a composite or average of the 
actual values of the adjacent layers. This effect is referred to as the shoulder 
effect. Even in the 0° case shown in the Figure 3A and 3B, where the tool is 
10 vertical and the formation is horizontal the measured data is quite different 
from the desired representation of resistivity. As the dip increases, the effect 
is increased. 

Much work has been done on methods and equipment for processing 
logging tool data or signals to produce an accurate representation of 

15 formation parameters. This data processing process is commonly called 
inversion. Inversion is usually carried out in some type of computer. In the 
prior art system of Figure 1 r a block labeled "computing module" may perform 
some type of inversion process. The methods currently available to perform 
this processing are iterative in nature. The standard iterative methods have 

20 the disadvantage of being computationally intensive. As a result, the 
inversion must normally be carried out at computing centers using relatively 
large computers, which can deliver results of the inversion in a reasonable 
amount of time, and normally cannot be performed in computers suitable for 
use at the well site. 

25 An alternative processing method is the deconvolution method. This 

method is very fast and can be implemented at the well site, for example in 
the computing module of Figure 1. However, this method is based on linear 
filter theory, which is an approximation that is not always accurate. In 
deviated boreholes, the nonlinearity of the tool response becomes manifest, 

30 making the problem hard for the deconvolution method to handle. The 
deconvolution methods do not generate actual representations of the 
formation parameters, so they cannot be properly called inversion methods. 
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Early attempts to solve the inversion of log data problem used the 
parametric inversion method. This method is an iterative method that uses a 
forward solver and criteria, such as the least square inversion, to determine 
the best fit for the parameters of a predefined formation, usually a model with 
5 a step profile. However, if the actual formation does not conform to the 
predefined model, the output parameters determined by this method can be 
very far from the actual parameters of the formation. This is a consequence 
of the ill posed nature of the inversion problem which makes it highly non- 
trivial. 

10 A more current method for inversion of resistivity log data is the 

Maximum Entropy Method, MEM. In this iterative method, a test or proposed 
formation model is modified to maximize the entropy functional, which 
depends on the parameters of the formation. This method does not use a 
predefined formation and produces solutions of better quality. It is more 

15 efficient than the parametric approaches, but is still computationally intensive. 
It can be applied to any type of tool for which a forward solver is available. An 
example of the MEM method is disclosed in US Patent 5,210,691 entitled 
"Method and Apparatus for Producing a More Accurate Resistivity Log from 
Data Recorded by an Induction Sonde in a Borehole." 

20 In general, all of the iterative inversion schemes have essentially two 

parts. The first part is a forward solver that generates a synthetic log from a 
synthetic test formation which is a reasonable representation of a real 
formation. The test formation is an assumed generally square-step plot of a 
formation parameter, e.g. resistivity, versus depth, like the plots 10 and 20 in 

25 Figures 3A and 3B. The forward solver simulates the response of a selected 
logging tool to the test formation to generate the synthetic log. If the logging 
tool has multiple transmitter receiver sets or arrays, as illustrated in Figure 1 , 
a separate forward solution is needed for each set, since each set responds 
differently. The second part of the iterative method is a criterion to modify the 

30 test formation. The criterion is based on the difference between the synthetic 
log corresponding to the test formation and the real log data measured by the 
tool. After the test formation has been modified, a new synthetic log is 
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generated by the forward solver. This process is repeated iteratively until the 
difference between the synthetic log and the real log is less than a predefined 
tolerance. The output of the inversion algorithm is the parameters of the final 
test formation. These parameters are plotted versus depth to produce the 
5 desired log. It is the iterative nature of these methods which makes them 
computationally intensive. 

SUMMARY OF THE INVENTION 

The present invention overcomes the foregoing and other problems 
with a method for determining a formation profile surrounding a well bore 

10 wherein upon receipt of field log data for the formation surrounding the well 
bore a Jacobian matrix is generated responsive to the received information. A 
log response is calculated responsive to the determined Jacobian matrix and 
a determination is made of whether the calculated log response converges 
with the received field log data. If the log response does not converge with 

15 the received field log data, a quasi-Newton update is performed on the 
Jacobian matrix and a new log response determined. A new convergence 
comparison may then be performed. Once the calculated log response 
converges with the received field log data a formation profile based upon the 
log response is output. 

20 BRIEF DESCRIPTION OF THE DRAWINGS 

A more complete understanding of the method and apparatus of the 
present invention may be obtained by reference to the following Detail 
Description when take in conjunction with the accompanying Drawings 
wherein: 

25 FIGURE 1 illustrates an induction tool; 

FIGURE 2 illustrates a well bore; 

FIGURE 3a through 3b illustrate computer simulations of logs taken by 
a well sonde; 
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FIGURE 4 is a flow diagram illustrating a first embodiment for 
determining parameters of an earth formation surrounding a well bore 
according to the present invention; and 

FIGURE 5 is a flow diagram illustrating the determination of initial 
5 Jacobian matrix using a sliding window vector. 

DETAILED DESCRIPTION OF THE EXEMPLARY EMBODIMENTS OF THE 
INVENTION 

Generally, processing of log data from a well bore uses an iterative 

10 inversion scheme having essentially two parts. The first part contains a 
forward solver that generates a synthetic log from given test information. The 
second part contains criteria to modify the test formation. The criteria is 
based upon the differences between a synthetic log, corresponding to the test 
formulation, and the real log measured by the well tool. After the test 

15 formation has been modified, a new synthetic log is generated by the forward 
solver. This process is repeated iteratively until the difference between the 
synthetic log and real log is a less than a predefined tolerance. The output 
and the inversion algorithm are the parameters of the final test formation. It 
should be pointed out, however, that the repeated computation of the forward 

20 model in each iteration makes these methods computationally intensive and 
require a great deal of processing time. 

Referring now to FIGURE 4, there is illustrated a flow diagram 
describing a first embodiment of the present invention for more efficiently 
determining parameters of earth formation surrounding a well bore. The field 

25 log data is received at step 30 and initially processed to determine at step 35 
a conductivity or resistivity vector for the log data. An initial Jacobian matrix is 
generated at step 40 using a sliding window as will be more fully described 
with respect to FIGURE 5. The Jacobian Matrix is then used with the 
maximum flatness constraint to solve for a new set of formation resistivity 

30 vector at Step 45. A new log response is calculated using the new resistivity 
vector at Step 50. The difference between the recalculated log response and 
the recorded field log is determined at step 55. Inquiry step 60 determines if 
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10 



15 



20 



25 



the provided field log data converges with the newly calculated response, and 
if so, the iterations are ended and a calculated log response is provided as the 
output formation profile at step 65. Otherwise, a quasi-Newton update of the 
Jacobian matrix is performed at step 70 utilizing the newly generated log 
response and a new inversion is performed at step 45. This process is more 
fully described below. 

1. PROBLEM PARAMETRIZATION 

Since an induction tool directly measures the formation conductivity, 
this parameter is used throughout the inversion process instead of resistivity. 
At the end of the inversion process, it is converted to resistivity. Normally one 
has no a priori knowledge of bed boundaries. This lack of knowledge is 
compensated by discretizing the unknown formation into thin beds of equal 
thickness. This thickness can selectively take the values of \ where h b can 

be Vi ft., 1 ft. or any integer multiple of the distance between consecutive 
logging points when the relative dip between borehole and formation is zero. 
Otherwise, the beds will have a thickness equal to }\ /cos(<9) where e is the 
relative dip angle. If there are a total of N such beds, the problem of inversion 
becomes one of finding a specific conductivity vector: 



such that the predicted log response based on this profile best, matches the 
actual field log. 

In order to keep all conductivity values positive, a traditional 
transformation is introduced by replacing the conductivity values with their 
logarithmic values during the inversion process: 
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and <y° =[a' l ffj l ...aj] is some initial conductivity profile which is typically 

chosen to be a uniform profile. Similarly, if there are a total of M logging data 
points in the field log, this field data is combined into a column vector of the 
5 form: 

And the inverted conductivity profile is written as: 

10 

The measured apparent conductivity of a multi-coil induction tool at a 
given logging depth z in a formation with conductivity profile a can be written 
15 as: 

a a (z,a) = ~Z ZVn M = -iHr rt (z,q) = a Q (z,q) (5) 

where 

20 

V rt (z t <r) = voltage generated on the r-th receiver due to the f-th 
transmitter. 

z = logging depth; 

T = total number of transmitters in the subarray; 
25 R = total number of receivers in the subarray; 

K = tool constant for the subarray. 

Since Eq. 5 holds at every logging position z = z m {m = 1,2,.,.,/W), one can 
construct an apparent conductivity vector over the whole logging interval as: 

11 
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°a (?) = \Po (*l > <ll <T(Zl > <1 1 • >°a > ?)f " ' ( 6 ) 

Because cr fl (^) is a nonlinear function of the formation conductivity vector a 
and hence q , just as in any nonlinear problem, one needs to first linearize the 
5 response function. To this end, a a {q) is expanded into a Taylor series about 
a conductivity profile q as: 

<r a {q + x) = ° a {<lhA( q )x + 0(x 2 ) (7) 

10 where x is a perturbation vector representing a small change in the 
transformed conductivity vector. A is the M x N sensitivity or Jacobian matrix 
whose elements are given by: 

A mn = da i Zm ' q) > m = l,2,...M * = l,2,:..,Ar (8) 

15 

When the Euclidean norm of = Vj^x is small, the higher order terms can 

be neglected. Eq.7 then represents an approximate linear dependence of the 
measured apparent conductivity on the transformed conductivity vector q. 
Since analytical expressions for the derivatives in Eq.8 are not available, 
20 finite-difference approximations are used. For example, the following forward- 
difference expression can be used: 

A m „-- °° {Z - q) ~ a ° {Z - q \ m = l,2,...,M,n = l2,...,N, (9) 
? = [*. •••> ?/>-, 9w] T . 

25 = <?2. • - ;<1nJ , 

where 8 is a predetermined small increment in the transformed conductivity. 
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2. THE MAXIMUM FLATNESS INVERSION ALGORITHM 

The linearized Eq.7 can be rewritten as: 
5 Ax = b (10) 

where the direct error vector or direct misfit vector is defined as: 

10 

Since there are typically more logging data points than the number of 
unknown conductivity values, Eq.10 is a standard over-determined least- 
squares (LS) problem, which can be solved readily to give: 

15 x = [A T Aj ] A r b. (12) 

To reduce overshoot in the least-squares solution, usually a damping term X 
is introduced into Eq.10 to form the damped least-squares (DLS) problem. Its 
solution can be written as: 



20 



x = [a t A + A 2 i]~ 1 A T b (13) 



where / is the identity matrix. This solution of the perturbation vector x is 
added to the initial conductivity profile q to yield the transformed conductivity 
25 profile of q + x. The resultant conductivity profile is then obtained via 
transformation (Eq.2) as: 

^ + x B =^\» = l,2,...^ (14) 
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Using k as iteration index such that q k denotes the transformed conductivity 

profile vector and x k the perturbation vector at the /c-th iteration, one can 
rewrite Eq.14 into the following iterative formula: 

5 g w (15) 

In order to stabilize the least-squares solution, some sort of 
regularization process needs to be introduced to the aforementioned 
procedure. Among the various kinds of regularization methods, the method of 

10 maximum entropy has enjoyed wide success in various branches of science, 
from the sharpening of radio astronomy data to pattern recognition, to x-ray 
crystallography. The maximum entropy principle was originally proposed as a 
general inference procedure by Jaynes based on Shannon and Weaver's 
pioneering work on information theory. 

15 For discrete conductivity profile such as the case of earth formations, 

the entropy functional is written as: 



25 



S{cr} = -Z^ln -Z*- , C = 2>„ (16) 



N 

I 



20 The negative sign means that this entropy is the so-called negentropy since it 
is the negative measure of the configurational information in the conductivity 

profile. The necessary condition for entropy to be the maximum — = 0 

da 

leads to: 

a n =a 0 ,w = l,2 f ...,tf (17) 



Namely, the entropy functional reaches its maximum when there is absolutely 
no feature in the conductivity profile, i.e., a homogeneous profile will have the 
highest entropy while the introduction of any feature into the profile inevitably 
reduces the entropy. Thus, by constraining the least-squares problem with 



14 

DALLAS2 882154V1 52298-0001 6USPT 



Patent Application 
52298-0001 6USPT 
2001-IP-5416 



the requirement that entropy be the maximum, one ensures that the inverted 
conductivity profile thus obtained will be the flattest or smoothest possible 
profile that is consistent with the original field data. 

With the transformation in Eq.2, the maximum entropy condition (Eq.17) 
5 is expressed in terms of the transformed conductivity vector q as: 

f = 0 (18) 

Since q k+l = q k + x k , setting q M =0 results in the following maximum entropy 
10 constraint: 

x = -q (19) 

The utilization of a Lagrange multiplier a enables one to combine the 
15 damped least-squares problem Eq.10 and the maximum entropy constraint 
(Eq.21 ) together to form a joint problem whose solution can be expressed as: 

x = [A T A + a 2 I + A 2 1 J' [A T b -a 2 q] (20) 

20 In essence, the starting point of the maximum flatness algorithm is the 
assignment of q n = 0 (w = 1,2,. ..,#). Namely, a homogeneous conductivity 
profile (see Eq.2) is assumed. This homogeneous conductivity profile has a 
maximum entropy or maximum flatness according to Eq.17. In subsequent 
iterations, the profile is modified so as to minimize the misfit W = |^-cr*|| 

25 between the field data a A and the predicted data a k . 

3. INVERSION IN TERMS OF RELATIVE MISFIT VECTOR 

The misfit vector defined in Eq.11 is directly proportional to conductivity 
30 and leads to inversion results that are biased toward highly conductive zones. 

15 

DALLAS2 882154V1 52298-0001 6USPT 



Patent Application 
52298-0001 6USPT 
2001-IP-5416 



Alternatively, one can introduce a relative error vector or relative misfit vector 
b r with the subscript Y' to distinguish it from the direct error vector b defined 
in Eq.1L Elements of the relative error vector are those of the direct error 
vector normalized by apparent conductivity: 

5 

^ g ^(0-cr (r., g ) a ^ m = la _ M (21) 



It is straightforward to show that elements of the new Jacobian matrix for this 
relative misfit vector are given by: 

10 

A rmn = / da ^ q \ m = \2^M\ n = l2,...,N (22) 



The effect of the relative misfit vector b r on the inversion is opposite to 
that of the previously defined direct misfit vector b in that the former is biased 

15 toward regions of higher resistivity while the latter is biased toward more 
conductive regions. Which algorithm should be used for a given field log 
should depend on the quality of the field log. If the field log data is "clean" in 
the sense that noise is known to be low, then the relative error algorithm 
usually yields better results. On the other hand, since log data in highly 

20 resistive beds are more susceptible to noise, for noisy data, the direct error 
algorithm is more robust against noise. 



4. QUASI-NEWTON UPDATE OF THE JACOBIAN MATRIX 



25 Because of the requirement of the Jacobian matrix, all algorithms 

discussed so far suffer the same problem of high computational cost. 
Specifically, given a total of M logging points, a total number of N discretized 
beds in the formation, and the total of /w. iterations, the total number of 
necessary calls to the forward routine is given by 
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n r = (MxN+M)k^ =M(N + \)k mas *MxNxk„. (23) 

Since the most time-consuming portion of the iteration is the calculation of the 
Jacobian, if one can reduce this part of the inversion, the whole algorithm can 
5 be accelerated significantly. A quasi-Newton updating technique is applied to 
the Jacobian matrix to reduce the inversion time. 

Equation 7 can be rewritten as: 

A M x = y k (24) 

1 0 where 

x*=fl* +, -f* (25) 



Instead of using the Jacobian matrix A, a matrix update B M = B k +U k is 
1 5 used that satisfies the quasi-Newton condition: 

B M x k =y k (26) 

Broyden's rank-one update formula is used to obtain the update in one 
20 embodiment of the current invention: 

rx*T 

= B k L J (27) 

M M 

The advantage of the quasi-Newton update algorithm over the direct 
25 Jacobian algorithm is obvious. In the quasi-Newton update algorithm, the 
Jacobian matrix is calculated only once at the very beginning of inversion 
process. In subsequent iterations, the Jacobian is updated with little 
computational expense since the data needed for this update is already 
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available after the forward modeling that is needed to calculate the misfit 
vector. Thus, the total number of calls to the forward modeling routine is: 

iQjt.-MxKn. (28) 

5 

Compared with Eq.23, the algorithm with quasi-Newton updates reduces the 
number of calls to the forward routine by a factor that is equal to the total 
number of logging data points. Thus, with thousands of logging points for any 
given section of field log data, quasi-Newton update result in huge savings in 

1 0 CPU time which makes fast inversion possible. 

It should be pointed out that the so-called rank-one quasi-Newton 
update is mentioned in the above. It is obvious that other types of 
quasi-Newton update methods can be used instead. These methods include, 
for example, the symmetric rank-one update, or the family of rank-two 

15 updates, such as the Davidon-Fletcher-Powell (DFP) update, 
Bryoden-Fletcher-Goldfarb-Shanno (BFGS) update. 

5. SPEED UP OF THE CALCULATION FOR THE INITIAL JACOBIAN 
MATRIX 

20 

Referring now to FIGURE 5, although the quasi-Newton update speeds 
up the inversion in subsequent iterations by avoiding direct calculation of the 
Jacobian matrix, before the start of the iteration procedure, the Jacobian 
matrix still needs to be initially populated. Fortunately, with the application of 

25 the maximum flatness algorithm to the field log data, the initial conductivity 
profile is a homogeneous one. This enables one to greatly reduce the 
number of calculations needed. Specifically, Except at the two beds at the 
two end points of the formation of interest where only a 2-bed formation is 
needed, the response of the induction tool to a simple 3-bed formation is 

30 needed to populate the whole Jacobian matrix. The 3-bed formation consists 
of a center bed having a thickness equal to the parameterized bed thickness 
(discussed in section 2) and two semi-infinite shoulder beds. With the 
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conductivity of each bed set at step 90 to the average of the logging data 
apparent conductivity, the response of the induction tool is calculated at step 
95. This necessitates only a single call to the forward routine. Then the 
conductivity of the center bed is incrementally changed at step 100 and the 

5 response of the tool is calculated at step 105. Forward difference (Eq.7) is 
then used to calculate the derivatives at step 110, which forms a single 
column vector of the Jacobian matrix. By sliding the 3-bed formation and this 
single column vector across the whole formation of interest at and keeping 
track of all the indices meticulously, the whole Jacobian matrix is populated at 

10 step 115. Compared to the brute force forward difference method, this 
technique reduces the CPU time required for the initial Jacobian matrix by a 
factor of M/ 3 where M is the total number of logging data points. 

The previous description is of a preferred embodiment for implementing 
the invention, and the scope of the invention should not necessarily be limited 

15 by this description. The scope of the present invention is instead defined by 
the following claims. 
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